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We present the results of Direct Numerical Simulations (DNS) of turbulent flows seeded with 
millions of passive inertial particles. The maximum Reynolds number is Re\ ~ 200. We consider 
particles much heavier than the carrier flow in the limit when the Stokes drag force dominates their 
dynamical evolution. We discuss both the transient and the stationary regimes. In the transient 
regime, we study the growth of inhomogeneities in the particle spatial distribution driven by the 
preferential concentration out of intense vortex filaments. In the stationary regime, we study the 
acceleration fluctuations as a function of the Stokes number in the range St £ [0.16 : 3.3]. We also 
compare our results with those of pure fluid tracers (St = 0) and we find a critical behavior of 
inertia for small Stokes values. Starting from the pure monodisperse statistics we also characterize 
' polydisperse suspensions with a given mean Stokes, St. 

I. INTRODUCTION 



Suspensions of dust, impurities, droplets, bubbles, and other finite-size particles advected by incompressible turbu- 
lent flows are commonly encountered in many natural phenomena and industrial processes. These finite-size particles, 
whose density may differ from that of the underlying fluid, cannot be modeled as point-like tracers because of their 
, inertia. Often, they are characterized by the presence of strong inhomogeneities in their spatial distribution. Indeed, 
light (heavy) particles tend to concentrate into specific regions of the flow characterized by high (low) values of the 
vorticity. Such a phenomenon is dubbed 'preferential concentration' The inhomogeneities then appearing in the 
particle spatial distribution are important because they affect the probability to find close particle pairs and thus 
influence their possibility to collide, or to have biological and chemical interactions. Examples showing the importance 
of this phenomenon are the formation of rain drops by coalescence in warm clouds 0j S S 1 or the coexistence 
between several species of plankton in the hydrosphere 0, 0^, Engineering applications encompass optimization of 
spray combustion in diesel engines and in rocket propellers!! ■ Inertial particles are also important for the problem 
of dispersion of dust, chemicals or aerosols in the atmosphere [lAlTl). 
, Recently, much effort have been devoted to the study of phe nomena related to preferential concentration of inertial 
particles in turbulent flows by means of both theoretical [12], [H El [13 and computational M E M E3 
rN | approaches. Progresses in the characterization of the statistical features of particle clusters have been achieved by 
studying inertial particles evolving in laminar stochastic flows 0, 152, IU l24l 1251 and two dimensional turbulent flows 
*26|. Experimental results are reviewed in Ref^Q. Recently, experiments |27| have been also realized borrowing 
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techniques from Lagrangian turbulence 28, 29, 30l l3lll32| . 

Considerable less attention has been paid to other aspects of inertial particle dynamics, and in particular to the 
statistical properties of their acceleration in turbulent flows. This is an important issue relevant, for example, to 
the development of stochastic models of particle motion |3Jj. It is also interesting to contrast the statistics of the 
acceleration of inertial particles with that of fluid tracers which has been extensively studied by means of experiments 
[Hii H3, Him, numerical simulations pi l3l l3^[3l Eoj and theoretical approaches [Hi El El This is even 
more interesting on the light of recent experiments |27j designed to study inertial particles in turbulent flows by means 
of Lagrangian measurements. Concerning fluid tracers, it has been found that the acceleration probability density 
function (pdf ) displays the typical fat tails of highly intermittent signals j^, [j^] • The tails of the pdf have been 
associated to trapping events into vortex filaments [23, Moreover, a model based on the multifractal description 
of turbulence has been proved able to describe the pdf of the acceleration tracers |3(| . Similar questions are still open 
for inertial particles. Because of the centrifugal force acting on such particles, inertia is expected to considerably 
change the statistics of trapping events into vortex filaments. It is important to quantify this effect by looking at the 
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modification of the pdf tails as a function of the degree of inertia. Moreover, it is not known how to generalize the 
multifractal approach used for passive tracers to the case of passive particles with inertia. 

In this paper we focus on the statistical properties of particles much heavier than the carrier fluid. We are interested 
in investigating, by means of high- resolution Direct Numerical Simulations (DNS), the combined effects of preferential 
concentration and of inertia in determining the statistics of acceleration. With this aim, we review and extend the 
results presented in Ref. -"- n particular, we perform a systematic study as a function of the Stokes number. Here 
Stokes is defined as St = Ts/t^, being r s the typical response time of the particle and t v the Kolmogorov time of the 
flow. 

Most results concern monodisperse suspensions, i.e. set of particles with a unique Stokes number. Since laboratory 
and natural suspensions are typically polydisperse, we shall extend our analysis to this case also. 
Understanding the relevant time scales to reach a stationary regime is important both for experimental studies - often 
limited in running time - and for a dynamical description of preferential concentration. We present an analysis of the 
transient regime, i.e. the time window necessary for the particles to reach a statistically steady regime, starting from 
their initial configuration. 

The paper is organized as follows. In Sec. ^ we recall the equations of motion for the case of inertial particles 
much heavier than the surrounding fluid. A detailed description of the numerical set up is given in Sec. lIIII In Sec. II VI 
the validity of the model equations for our numerical experiments is considered. Moreover, we discuss the statistical 
features of preferential concentration both during the transient and in the statistically stationary regime. The most 
relevant features of particle acceleration statistics for monodisperse suspensions are discussed in Sec. [V] Section IVTl 
is devoted to extend the previous analysis to polydisperse cases. Conclusions and perspectives are drawn in Sec. VII 



II. HEAVY PARTICLE DYNAMICS 



The equations of motion of a small, rigid, spheric al p article immersed in an incompressible flow have been consis- 
tently derived from first principles in Refs. |44i l45l l4q under a certain number of assumptions. The main working 
hypothesis is to consider low particle concentration, so that it is possible to neglect collisions, particle-to-particlc 
hydrodynamic interactions and the feedback of particles on the carrier fluid. Another hypothesis is that the parti- 
cle radius r is much smaller than any active scale of the turbulent flow, namely that r <C 77, where rj denotes the 
Kolmogorov dissipative scale. Finally, the particle Reynolds number Re p — r\V — u\jv, where \V — u\ is the typical 
relative velocity of the particle with the flow and v is the fluid kinematic viscosity, has to be small, i.e. Re p <C 1. 
These hypotheses lead to describe the fluid flow surrounding the particle as a Stokes flow and to write a close dif- 
ferential equation for the motion of the particle |47| . The dynamics of a single particle then depends on only two 
dimension-less parameters. The first is = 3pf/(pf + 2p p ), accounting for the mass density ratio between the particle, 
p p , and the fluid, pf. The second control parameter is the Stokes number St = t s /t v , where the particle response 
time is r s = r 2 /{3[3v). Further, neglecting gravity, and assuming that the particle is much heavier than the fluid, the 
particle velocity is to leading order a solution to 0, : 

dV Rn V(t)-u(X(t),t) 

— = (JD t u(X,t) . (1) 

at t s 

where X denotes the particle position, u(x,t) the fluid velocity field, and D t u(X,t) — dtu + u ■ Vu is the fluid 
acceleration at the particle position. As long as D t u(X,t) does not become too large (see Sect. IIV All . the first term 
in the right-hand side of can be neglected when (i <C 1, so that the equations can be further simplified to 

¥ = ~-{V{t)-u{X{t),t)) . (2) 
at t s 

It is interesting to compare the evolution equations for particles (|2J) with that for the motion of fluid tracers 

dx(t) 



dt 



= u(x(t),t), (3) 



which corresponds to the limit r s = of (J2J. Tracers may be seen as St = particles. 
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III. DETAILS ON THE DNS AND THE DATA SET 

The carrier fluid is evolved according to the incompressible Navier-Stokes equations 

du „ Vp . „ 

— +u-Wu = + vAu + f, (4 

dt Pf 

where p is the pressure field and / is the external energy source, (/ • u) = e. These are solved on a cubic grid 
with periodic boundary conditions. Energy is injected by keeping constant the spectral content of the two smallest 
wavenumber shells |48|. The viscosity is chosen so to have a well resolved dissipative range, r\ « Ax where Aa: is 
the grid spacing. We use a fully dealiased pseudospectral algorithm with 2 nd order Adam-Bashforth time stepping. 
We performed three sets of runs at resolution N 3 with N = 128, 256, 512, the corresponding Reynolds numbers are 
Rx = 65, 105, 185 (see also Ref. 0)- 

Particle dynamics, Eqn. ©, is integrated in parallel with that of the fluid by means of a 2 nd order Adam-Bashforth 
time-stepping, in which the fluid velocity at particle position is estimated by linear interpolation. 

For each set of runs, the following procedure has been adopted. We integrate the Navier-Stokes equations until the 
flow reaches a statistically steady state. Once a turbulent configuration for the flow is obtained, millions of particles 
and tracers are seeded in the flow. Their initial positions are uniformly distributed in the volume, and velocities are 
chosen equal to the local fluid velocity. Equations J3J) for particles and © for tracers are then advanced in parallel 
with those for the fluid We consider particles with 15 different Stokes numbers ranging from 0.16 to 3.5 (see 
caption of Table [I] for the exact values), and a set of tracers, St = 0. For each set of particles with a given Stokes 
number, we save the position and the velocity of N t particles every dt = r, ; /10, with a maximum number of recorded 
trajectories of Nt = 5 • 10 5 for the highest resolution. Along these trajectories we also store the velocity of the carrier 
fluid. At a lower frequency, ~ 10t i( , we save the position and velocity of a larger number N p of particles (up to 7.5 • 10 6 
per St, at the highest resolution) together with the Eulerian velocity field. 

A summary of the physical parameters for the run at the highest resolution, the one discussed in this paper, is 
given in Tabled 



R\ Uims £ v V L T E 7V, 

185 1.4 0.94 0.00205 0.010 n 2.2 0.047 



Ttot Ttr Aa: N 3 N t N p Ntot 
12 4 0.012 512 3 5 ■ 10 5 7.5 ■ 10 6 12 ■ 10 7 



TABLE I: DNS parameters. Microscale Reynolds number R\, root-mean-square velocity u rms , energy dissipation e, viscosity 
v, Kolmogorov length scale r\ = (u 3 /e) 1 / 4 , integral scale L, large-scale eddy Eulerian turnover time Te = L/w rms , Kolmogorov 
timescale t v = [v/e) 1 ^ 2 , total integration time Ttot, duration of the transient regime Tt r , grid spacing Ax, resolution N 3 , 
number of trajectories of inertial particles for each Stokes Nt saved at frequency T n /10, number of particles N p per Stokes 
stored at frequency lOr^, total number of advected particles Ntot' Errors on all statistically fluctuating quantities are of the 
order of 10%. The Stokes number are St = (0; 0.16; 0.27; 0.37; 0.48; 0.59; 0.69; 0.80; 0.91; 1.01; 1.12; 1.34; 1.60; 2.03; 2.67; 3.31) 



IV. STATISTICAL ANALYSIS OF THE DATA SET 

First, we check if the approximate equations @ are consistent with the assumptions used to derive them, in the 
range of parameters of our DNS. Second, we study how particles, injected homogeneously in the domain, reach a 
statistically inhomogeneous steady state spatial distribution. Further, we characterize preferential concentration of 
particles as a function of the Stokes number in the stationary regime. 

A. Limit of validity for the particle equations of motion 

To validate our working hypotheses, we check that they are self-consistent with neglecting the term proportional to 
(3 in QJ. in particular, we want to see what are the constraints on the values of /3 and on the particle radius r that 
should be satisfied in order to fulfill the assumption (3\D t u\ <C \u — V\/t s . 
The first observation is that from the definition of St and from the dimensional estimation of the Kolmogorov scale 
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FIG. 1: Log-log plot of the pdf P(Q) vs Q = f3r s \dtu\/\u - V\, with (3 
St = 0.16, 0.59, 3.31 from the inner to the outer curve, respectively. 



0.001 and for three different Stokes numbers 



7/ = (^ 3 /e) 1 / 4 , one has rjr\ = y/3f3St , meaning that for St <S 1 there are no severe constraints on the value of f3. 
Differently, for large St, in order to have r <C one may need to use unphysical small values of /3. To see whether 
this is the case, let us consider a typical example of parameter values: the case of water in air (3 w 0.001, and St = 3 
(which is the largest value used in our simulations). In the worst case, we have r/rj < 0.3 which is still small enough to 
justify to neglect the Basset history and the added mass terms, and the Faxen corrections 0], leading to Eq. for 
the particle dynamics. We now show that it is self-consistent to further neglect the term (3D t u. Since in the DNS we 
do not have stored the fluid acceleration at particle positions for the whole trajectory, we assume D t u ~ dtu, i.e. we 
estimate the fluid acceleration at the particle position with the time derivative of the fluid velocity along the particle 
path. The validity of this approximation relies on the fact that the pdf of dtu that, for large St, are very close to 
those of D t u (not shown), meaning that strong acceleration events are well represented. 
Figure ^ shows, for three representative Stokes numbers, the distribution of the ratio 



Q = P 



\du/dt\ 
\{u-V)\/i 



(5) 



As one can see, the statistical weight of events Q > 0.1 is very low, meaning that with a very high probability the 
equations @ are a valid approximation. Of course, for density ratios smaller than f3 = 0.001 the approximation will 
work even better. If one considers for instance sand dust transported in the atmosphere at an altitude of roughly 
5km, (3 is reduced by a factor 1/10 and in this case our simulations do not catch any event where Q > 0.1. 

We thus conclude that in the range of (3 considered above, the fluid inertia term is also negligible and equations 
correctly describe the motion of small heavy particles in a turbulent flow at moderate Reynolds number. 



B. Transient regimes and formation of preferential concentrations 



It is well known that inertial particles tend to concentrate into specific regions of the flow. By performing a local 
analysis of the equations of motion it is possible to show that, for small Stokes numbers, particles heavier than the 
fluid tend to concentrate onto strain dominated regions of the flow and escape from rotational regions, while the 
opposite is observed for lighter particles |15l llo . 12 1| . Hence, particles do not sample uniformly the full velocity field 
and are sensitive to its geometry. 

Being interested in those properties of particle dynamics that are closely linked to the flow geometry at the particle 
locations (such as the distributions of velocity and of acceleration), one should clearly understand the bias induced by 
the correlation between particle positions and flow structures. Indeed, highly intermittent spots of the turbulent flow 
have a strong signature on the acceleration distribution of simple tracers through the phenomenon of vortex trapping 
[28l |35|. Moreover, at variance with tracers, whose dynamics is completely determined by that of the velocity field, 
inertial particles have their own dynamics. Therefore, for arbitrary initial data, a relaxation time is needed for them to 
reach a statistically stationary motion. Our simulations are started with a spatially uniform distribution of particles 
with a vanishing acceleration (velocity equal to that of the fluid). At sufficiently large times, the spatial distribution 
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FIG. 2: (Left) For = 0.16, particles positions on a two dimensional slice of the simulation box and, in grey scale, the 
amplitude of the vorticity field. (Right) The same for St = 2.03. 



of particles becomes strongly inhomogeneous and strongly correlates with the flow structure (see Fig.|21for two typical 
examples of the particle spatial distribution). As one can see, even for small values of the Stokes number, particles 
get less concentrated in those regions where the flow develops strong vortical structures, represented as dark grey 
spots on the figure. To have insights into the dynamics of heavy particles and to control the statistical analysis it is 
important to understand the time scales involved in the dynamical process of formation of such inhomogeneities. 

A naive inspection of the equation for the motion of a single particle suggests that after few r s , the particle 
velocity should relax to that of the fluid. However this relaxation time is clearly not enough to stabilize the statistical 
properties of the particle distribution, and in particular to form statistically stationary clusters of particles. Indeed, 
at least for large St, particles organizes in structures with length scales comparable with the largest ones present in 
the system (as shown in the right panel of Fig. [2J • 

To quantify the formation of inhomogeneities in the particle distribution we used two observables. The first one 
gives a rather global insight and is based on the coarse-grained mass distribution of particles. The second observable 
is local and is related to the particle distribution at scales of the order of the most violent events of the turbulent 
flow, namely the fraction of particles in rotation and stretching regions of the flow. 

We start the characterization of the inhomogeneous distributions of particles by measuring the degree of concentration 
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FIG. 3: (m 2 ) vs time, for particles with Stokes numbers St = (bottom solid line), St = 0.16 (empty boxes), St — 0.27 (filled 
boxes), St = 0.48 (empty circles), St = 0.90 (filled circles), St = 1.60 (empty triangles) and St — 3.31 (filled triangles). The 
two vertical arrows indicate the large scale eddy turnover time Te ~ 2. (dashed arrow) and the time chosen to be the end of 
the transient dynamics T tr ~ 4 (solid arrow). 



as obtained by coarse-graining the volume in small cubes of side A, and then looking at the distribution of the number 
rrij (or mass) of particles in the j-th cell. Figure |3] represents the behavior of the second-order spatial moment 
(to 2 ) = (A/2-7r) 3 ; m | °f the mass distribution measured on cubes of side A = 2Ax. For tracers, which remain 
uniformly distributed, (to 2 ) stays constant to the initial value, as given by a Poisson distribution. Differently, as soon 
as St > 0, (to 2 ) starts to increase in time until it stabilizes, fluctuating around a mean value. As shown in the figure 
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all curves reach this steady state only for (>2« i.e. for times of the order of the large scale eddy turnover time. 
In performing statistical analysis, for the sake of safety, we consider that the transient ends at T tr ss 2Te (second 
arrow in the figure). Let us notice however that, even in the statistically steady regime, t > T tr , fluctuations of (m 2 ) 
can be rather large and correlated to those of the total energy (see Fig. , making the definition of the transition 
between unsteady and steady regimes somehow difficult. At lower resolutions, we obtained quantitatively similar 
behaviors. We also checked that the time scale for equilibration does not depend sensibly on A (this remains true 
only if A is not too large). 

We now turn to the second observable based on the particle position conditioned on the flow local geometry. There 




FIG. 4: Time evolution of the total energy E = | J dr\u(r, t)\ 2 (continuous red line) and of the total energy dissipation, e 
(dotted blue line). 

are several possible ways to identify strain or rotation dominated structures in three dimensions (see Ref. |49( for a 
review) : here we characterize the geometry of the fluid velocity field u by looking at the eigenvalues of the strain 
matrix <Xjj = diUj. Generally, if the eigenvalues are all real, the point is said to be hyperbolic, while if they are all 
imaginary it is elliptic. In two-dimensions hyperbolic and elliptic points clearly identify strain and vortical structures, 
respectively. In 3d, pure elliptic points do not exist and the identification is less straightforward; the hyperbolic or 
non-hyperbolic nature of a point can be then identified by looking at the sign of the discriminant [50| : 

In deriving JfjJ), we omitted the term proportional to Tr[er] because of incompressibility. For A < the strain matrix 
has 3 real eigenvalues (strain regions), while for A > it has one real and 2 complex conjugate eigenvalues (rotational 
regions) . The number of particles N$t ( A > 0) in A > regions is represented in Fig. as a function of time and 
for different Stokes numbers. We chose here to normalize this quantity by the number iVs t= o(A > 0) of tracers in 
non-hyperbolic regions to smooth out the instantaneous fluctuations of hyperbolic and non-hyperbolic points in the 
flow. At t = this ratio is 1 for all values of St because particles are uniformly injected in the domain. As time goes 
on, it follows a very fast drop from the initial value. Unfortunately, the strong fluctuations, correlated to those of 
the total energy (Fig. 0J, do not permit to unambiguously clarify if the time required for Nst(A > 0)/ATst=o(A > 0) 
to stabilize is or not of the same order of that needed by (m 2 ). However, according to Fig. El it does seem that 
the correlation with the flow geometry settles on a time scale shorter than that needed to stabilize particles mass 
moments. It would be interesting to test this point, by studying the simultaneous evolution of particles uniformly 
injected, and particle already stabilized with the surrounding flow. This would help to disentangle spurious effects 
induced by the fluctuation of the global properties of the flow. 

We conclude this section by commenting the right panel of Fig. where we show the average number of particles 
normalized by the total number N p , (Ag t (A > 0)) /N p in the stationary regime. As previously observed |43j . it appears 
that preferential concentration is non-monotonic in St, but there is an "optimal" Stokes number for which the effect 
is maximal, here St ~ 0.55. This behavior can be understood as follows. To have particles concentrated into specific 
flow structures, their response time must be fast enough to allow them to follow the flow evolution. Clearly large St 
particles are unable to do that and tend to decorrelate from the flow geometry; though their large scale distribution 
may still be strongly inhomogeneous (as indeed observed, Fig.HJ. On the other hand, for St<^l particles behaviour 
is closer and closer to that of tracers, whose positions do not correlate with the flow geometry. As we will see in 
next section the balance between these two effects has a strong signature in the probability distribution of particle 
acceleration. 
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FIG. 5: Average fraction of particles in the non-hyperbolic regions Nst(A > 0)/iVst = o(A > 0) vs time, for Stokes numbers 
St = 0.16 (plus), St = 0.27 (crosses), St = 0.48 (stars), St = 0.90 (empty boxes), St = 1.60 (filled boxes) and St = 3.31 (empty 
circles) . The vertical arrows are as in 

Fig. El Right: average fraction of particles in the non-hyperbolic regions {Nst(A > 0))/N p 
vs St. The average has been performed considering only steady state data, i.e. for t > T tr . 



V. STATISTICAL PROPERTIES OF THE ACCELERATION 



One of the most striking features that characterize tracers is the very intermittent distribution of the acceleration 
HE HE EE HH : fluctuations as strong as 80 times the root mean square acceleration a rms reflect the tendency of 
tracers to be trapped into vortical structures HE IBH l5^ | . We would like to understand if this property is shared also 
by inertial particles and what are the main features of the acceleration fluctuations, at changing the importance of 
inertia. A natural way to proceed is to start from the equations of motion. In the statistically stationary state, the 
formal solution of Eq.J5J is: 

V(t) = - / e- {t - s V T °u{X{s),s)ds 

T s J-oo 

1 f* 

a(t) = — e- {t - s)/T °[u{X{t),t) -u(X(s),s)]ds . (7) 

It is clear that the kernel e~^~ s ^ Ta acts as a low-pass filter on the fluid velocity differences, suppressing frequencies 
larger than the inverse of the Stokes time t~ . The larger is the Stokes number St, the more important this effect 
should be. On the other hand, we have indications that the limit of vanishing Stokes number is singular, due to the 
clusterization of particles onto fractal sets p3 . l2fj| . This should have consequences on the acceleration statistics also. 
In Fig. the probability density functions for the particles acceleration are plotted for different Stokes numbers and 
compared with that of tracers. As we could expect, acceleration fluctuations become less and less intermittent as 
St increases. To disentangle different effects, it is useful to analyze the two limiting cases of small and large Stokes 
numbers, separately. 

At small St, i.e. r s <C T, q , the fluid velocity along the trajectory evolves smoothly in time and the acceleration 
reduces to a(t) ~ d t u(X (t) , t) , i.e. to the derivative of fluid velocity along the inertial particle trajectory. In the 
limit of vanishing St, this is indistinguishable from the fluid acceleration D t u(X(t),t) evaluated at particle positions 
and one expects that the particle acceleration essentially coincides with the fluid acceleration. There is however one 
major difference with the tracer case : inertial particles preferentially concentrate inside regions with low vorticity 
(see right panel of Fig. - In Fig- El we plot the results for a rms , as a function of St. Notice the singular effect of the 
inhomogeneous spatial distribution which results in a drastic reduction of a rms already for small St. This behavior 
is even more pronounced for the fourth moment of the acceleration |43| . To better understand the importance of 
preferential concentration, we have measured the tracers acceleration ((Dtu) 2 ) 1 ^ 2 conditioned on the spatial positions 
of the inertial particles: the result is also plotted in Fig. [7| The agreement of the two curves at small St confirms the 
validity of the previous argument. 

However, at increasing St, the curves start to deviate. While a rms monotonically decreases in St, the tracer 
acceleration conditioned on the particle positions has a minimum for St w 0.5 close to the maximum of preferential 
concentration, eventually recovering the value of the unconditioned tracers for larger St. For increasing St, we 
indeed know that inertial particles explore the small scale structures of the flow more and more homogeneously and 
a different mechanism is responsible for the reduction of a rms . 

At large St, i.e. r s 3> r n , we indeed expect the filtering effect to become more and more important. Again, we start 
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FIG. 6: Normalized acceleration pdf's for a subset of Stokes values (St = 0, 0.16, 0.37, 0.58, 1.01, 2.03, 3.31 from top to bottom) 
at R\ = 185. The darker line corresponds to tracers, St — 0. 




FIG. 7: Comparison between a rms (fill circles) as a function of Stokes with the acceleration of the fluid tracer conditioned to 
be on the particle position, ((Dtu) 2 ) 1 ^ 2 (empty boxes). The last curve (empty circles), approaching a rms for large St, is the 
one obtained from the filtered tracer trajectories, a^ ms . 



from tracers: for each solution, x(t) of Eq. Q, we define a new velocity, u F , obtained from filtering the fluid velocity 
along the particle trajectory over a time window of the size of the order of the Stokes time: 

1 r* 

u F (t) = - e- {t - s)/T °u(x(s),s)ds, (8) 

T s J-co 

then the filtered acceleration is given by a F — du F jdt. In figure the acceleration variance of particles a rms is 
compared with that obtained from tracers in the above manner, without any additional spatial conditioning. The 
curves corresponding to a rms and to af ms become closer and closer as St increases, supporting the conjecture that 
preferential concentration for St > 1 becomes less important. 

For intermediate St a non trivial interplay between the two above mechanisms takes place. A model, even qualitative, 
able to bridge the gap between the two limits in still unavailable. 

In the limit of small and large St the qualitative trend of the pdf's can be captured by the same arguments as 
those used for a rms . In Fig.|S]we compare the normalized pdf particle acceleration with those obtained by using the 
tracer acceleration measured on the particle position, D t u, and the filtered tracer trajectories, a F , respectively. For the 
smallest Stokes number St = 0.16, the overlap between the conditional tracer acceleration and the particle acceleration 
is almost perfect, while the filtered one is very close to the unconditioned tracer acceleration pdf. The opposite is 
observed for the largest Stokes number, St = 3.31. Notice that the (normalized) conditional tracer acceleration pdf 
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FIG. 8: In each panel the following four curves are shown: the normalized acceleration pdf for inertial particles (solid line), 
the acceleration pdf of fluid tracers Dtu conditioned to particle positions (empty boxes), the acceleration pdf of filtered fluid 
trajectories a F (filled boxes) and the unconditioned acceleration pdf for fluid tracers (dotted line). Note that the filtered 
acceleration a F pdf, which is for small St practically coincident with the (unfiltered and unconditioned) fluid acceleration pdf, 
becomes closer and closer to the particle acceleration pdf as St increases. The opposite is observed for the tracer acceleration 
conditioned to particle positions. 



is very close to the unconditioned one, meaning that preferential concentration is no longer effective on the small 
scale statistics for large St. For two intermediate St values, we observe the transition from preferential concentration 
dominated regime to the filtering dominated one. 

We conclude this section by going back to the instantaneous distribution of particles for small and large St (see 
Fig. El- Both pictures present an inhomogeneous spatial distribution for the inertial particles but for different reasons. 
At small St (left panel), inhomogeneity is mostly related to the fact that particles escape from large vorticity regions. 
For such small response times t s , particles are able to follow the flow variations but the very small scales properties of 
the fluid are sampled in a biased way. This explains the remarkable agreement between particle acceleration statistics 
and fluid tracers conditioned on particle positions. At larger Stokes number, the inhomogeneity is less correlated 
to the fine structures of the flow. Indeed we observe void regions whose typical size can be much larger than the 
Kolmogorov length scale and comparable to the inertial range scales. In spite of the clear inhomogeneous distribution 
of particles, the small scales, r <C rj, are sampled more homogeneously than at smaller St. This is coherent with the 
fact that the conditional tracer acceleration statistics is very close to the unconditioned one. 



VI. ACCELERATION OF POLYDISPERSE SUSPENSIONS 



Suspensions generally involves particles distributed over different sizes. This is true not only for phenomena observed 
in nature, as for aerosol distribution in the atmosphere, but also for laboratory experiments. In some cases, particles 
can also considerably change their size with time, as it happens for water droplets in warm clouds so varying their 
Stokes number. Direct numerical simulations allow to study inertial particles in the idealized case of a given Stokes 
number, as we have done previously, but also permit to test phenomenological approaches for polydisperse flows. In 
this section, we study what happens to the acceleration pdf when particles have a non trivial size distribution, or 
equivalently when the Stokes number of a given suspension is a random variable with a distribution Q(St). This 
is particularly important to have clean and possibly unambiguous comparison of numerical simulations results and 
laboratory experiments. 

To model the acceleration pdf V P oiy(a) in polydisperse flows, we first need a functional form for the pdf, V(a\St), 
of particles with a given Stokes number St. Second we make a convolution of the latter distribution with Q(St). 
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FIG. 9: Comparison of the acceleration pdf from DNS data (symbols) and the fit (solid line) of eq.® for two monodisperse 
cases: for St = 0.27 (outer curves) and St — 1.01 (inner curves). 
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FIG. 10: Comparison between the acceleration pdf for a polydisperse case (filled squares) of mean Stokes value St, and a 
monodisperse case (solid line) at the Stokes number St — St. In the inset, the Poisson distribution Q(St) of mean Stokes value 
St. Left panel refers to case with the Stokes number St = 0.27, while right panel refers to St = 1.01. 



Building up a phenomcnological model for the acceleration pdf of monodisperse suspension valid at different Stokes 
number is not an easy goal. For tracers, the multifractal formalism, without free parameter, is able to predict an 
acceleration pdf which is in excellent agreement with DNS data j3|| . However any simple extension of such approach 
to heavy particles does not seem to work, because of the almost singular nature of the small St limit and of the 
complex interplay between preferential concentration and time filtering. In the lack of such a model, the first step is 
done by fitting the acceleration pdf at changing the Stokes number. At small values of the acceleration, i.e. in the 
pdf core, we expect to have a Gaussian-like profile V(a\St) ~ exp(— a 2 ); while for large fluctuations a 3> a rms , we 
should recover a distribution with stretched exponential tails, V(a\St) cx exp(— a 7 ). A simple functional form for the 
pdf, V(a\St) of the normalized acceleration y — a/a rms (Si) satisfying these two requirements is: 

V{y\St) = a rms V(ya rms \St) = exp (- - Cl — + cA. (9) 

Here, to have a simple fitting procedure, we have used four free parameters Ci(St) with i = [1, ... ,4], although two 
of them may be fixed by the requirement of having the pdf normalized to unit area and unit variance. Of course, the 
most important parameter is C2(St) which fixes the Stokes number dependency of the far tail exponent 7 = 2 — c-i- 

The acceleration pdf measured in our DNS has been fitted with expression © for all available St. The values for 
the four free parameters are those reported in Table IVT1 They all gently, but non monotonically, vary with the Stokes 
number. In particular, the exponent C2, which controls the fat tails, goes from about ~ 1.6, for the smallest Stokes, 
to ~ 1.2, for the largest available one. This is in agreement with the tendency towards a less intermittent distribution 
for particles with a larger inertia. From Fig. EI we can see that the fitting formula ijjjj reproduces well the acceleration 
distribution (here shown for St — 0.16 and St = 1.01 only). 

To assess the variations for a polydisperse suspension, we make the convolution between the acceleration pdf condi- 
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St 





0.16 


0.27 


0.37 


0.48 


0.59 


0.69 


0.80 


St 


0.91 


1.01 


1.12 


1.34 


1.60 


2.03 


2.67 


3.31 


Cl 


2.68 


2.17 


2.24 


2.31 


1.67 


1.64 


1.83 


1.66 


Cl 


0.83 


1.64 


1.15 


1.08 


1.18 


1.95 


0.89 


0.76 


C2 


1.55 


1.43 


1.38 


1.35 


1.36 


1.33 


1.29 


1.28 


C2 


1.40 


1.25 


1.30 


1.28 


1.24 


1.11 


1.23 


1.27 


C3 


0.79 


0.75 


0.83 


0.91 


0.62 


0.63 


0.76 


0.68 


C3 


0.24 


0.69 


0.42 


0.39 


0.46 


1.00 


0.30 


0.21 


C4 


-0.41 


-0.53 


-0.53 


-0.53 


-0.65 


-0.67 


-0.64 


-0.68 


c 4 


-0.98 


-0.69 


-0.82 


-0.85 


-0.81 


-0.67 


-0.90 


-0.94 



TABLE II: The values of the fitting parameters d at changing the Stokes number. Notice that only two out of the four 
parameters are indeed independent because the pdf must be normalized to have unit area and unit variance. The fitted 
functional form is not meant to have a solid phenomenological base. It is used, as simple as it is, to estimate the importance 
of polydisperse suspensions with respect to monodysperse ones. 



tioned on the Stokes number V(a\St) and the Stokes number distribution Q(St) : 

V po i(a)= J dStV(a\St)Q(St). (10) 

A possible choice for Q(St) is a Poisson-like distribution: for any available Stokes number St, we generate a Poisson 
pdf with mean value St — St and we insert it as test curve to evaluate the convolution (|10ll . In Fig. ^|we compare the 
polydisperse pdf Ijl0|l obtained with a Poissonian distribution, Q(St), with a given St and a monodisperse distribution 
with St = St. We do it for two different average Stokes number, namely St ~ 0.27, 1.01. 

As expected, the effect of a random distribution for the Stokes number is much more important for the tails than 
for the cores. We recall also that, in the case of monodisperse suspensions the pdf tails are the most sensitive to St 
variations. Hence, any polydisperse suspensions characterized by a particle size distribution skewed towards small 
values of the Stoeks numeber should develop more intermittent tails than those expected on the basis of their mean 
Stokes value. 



VII. CONCLUSIONS AND PERSPECTIVES 



Direct numerical simulations, even at moderate Reynolds number, represent a valuable tool for the study of the 
Lagrangian motion of heavy particles. Here we have described the specific set up of our numerical study, analyzed 
the importance of preferential concentration and particle clustering as a function of the Stokes number. Moreover, 
we have explored the transient dynamics which precedes the settling of the particle motion onto a statistically steady 
state in terms of the particle distribution. This may be particularly important for real observations and laboratory 
experiments, where long time records of particle motion are hardly attainable. We have also reviewed and extended 
some recent results about the acceleration pdf of inertial particles, at varying the importance of inertia. The main 
conclusions are (i) preferential concentration plays an almost singular role at small Stokes, since even a small inertia 
is sufficient to expel particles from vortical regions where the strongest acceleration fluctuations are experienced; (ii) 
for small Stokes, a good quantitative agreement between the inertial particle acceleration and the conditioned fluid 
tracer acceleration is obtained; (iii) at large Stokes, the main effects is filtering of the velocity induced by the response 
Stokes times. 

The formulation of a phenomenological model able to describe the inertial particle acceleration as a function of both 
Stokes and Reynolds numbers is still beyond our reach. An easy fitting procedure has however been proposed for 
polydisperse suspensions, the most common experimental situation. For polydisperse cases, we show that the simple 
approach of characterizing the suspension in terms of its average Stokes number without taking into account the 
Stokes number distribution can give a systematic bias, particularly if the particle distribution is skewed towards small 
Stokes values. It would be extremely important to test the proposed approach with experimental data. Finally it 
would be interesting to numerically investigate St ^ 1 particles to get a better understanding of the almost singular 
behavior that has been observed for the acceleration. 
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